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MINIMUM WEIGHT DESIGN OF STRUCTURES 
VIA OPTIMALITY CRITERIA 

INTRODUCTION 


Every design task, whether applied to a structure or some other man- 
made object, is either directly or indirectly governed by certain optimality 
criteria. The designer is seldom required to create a product that will only 
serve its primary function — he is invariably expected to meet additional de- 
sign objectives, such as the lowest possible weight or cost, minimal main- 
tenance, maximum reliability, pleasing appearance, etc. 

In structural mechanics, the term "optimal design" is commonly used 
in a much more restricted sense. It implies that the sole design objective is 
minimum structural weight . 1 Moreover, it is understood that the entire design 
process is carried out automatically on a digital computer. 

The basic idea of computer-automated, minimum weight design can be 
best comprehended by viewing the optimizing algorithm as an evolutionary 
advance over the conventional design method. Consider a typical problem in 
which the designer is given the layout of the structure, the loading, and the 
constraints on the behavior of the structure (e.g. , limits on stresses and 
displacements). He must then calculate the sizes of the members so that the 
total structural weight is minimized and no constraints are violated. 

The design procedure is an iterative, trial-and- error process, each 
iteration consisting of two steps: an analysis of the current design, followed 
by a redistribution of structural material. In the conventional design method, 
the analysis, which we presume is carried out by a finite element computer 
program, is used primarily to check the behavioral constraints and provides 
little or no guidance for the material redistribution cycle. Consequently, the 
redesign is a creative task, based on empirical rules and the intuition and 
experience of the designer. 

An obvious means of improving the technique is to extend the scope of 
the analysis cycle so that it will not only calculate the behavior of the current 


1. Optimal weight design methods can also be used for minimum cost design 
by replacing the unit weight of each structural element by its unit cost. 



design but also predict the changes in the behavior caused by material redis- 
tribution. As a specific example, consider a structure subjected to displace- 
ment constraints only. Let us assume that an analysis of the structure indi- 
cates that some of the displacements exceed the prescribed limits so that the 
designer must add to the structural stiffness by increasing the sizes of some 
of the members. The question is, which members and by how much? 

To obtain a rational answer, the designer must have some idea of how 
sensitive the displacements are to changes in member sizes. This could be 
accomplished, at least partially, by computing the displacement gradients in 
addition to the displacements themselves. By displacement gradients we mean 

3uV 3W. , where u. is the i^ generalized displacement and W. represents 

th 

the weight of the j member. The designer would then know which members 
are most effective in reducing the displacements that violate the constraints 

and by using the approximation Au. = ^ (3u./3W.) AW. , he could also 
estimate the magnitude of the required redesign. 

In practical problems such redesign data are too voluminous for the 
designer to handle. It is much more reasonable to process the data within the 
computer itself, i. e. , to let the computer take care of the redesign in addition 
to the analysis. We have now arrived at a fully automated, optimal design 
algorithm. 

It should be pointed out, however, that it seldom is feasible to obtain 
an efficient design without any human participation whatsoever. Structural 
layout, for example, is one branch of design that is unlikely to be automated 
in the near future, with the exception of some specific problems (e.g. , trusses). 
Even optimal weight distribution with fixed layout, which is the topic of this 
report, cannot be entirely automated. In many problems it is neither prac- 
tical to include all the constraints into the computer program (the subjective 
constraints defy formulation altogether) nor is it possible to arrive at an 
efficient choice of member types on the first attempt. All this requires con- 
tinuous monitoring of the design process by an experienced engineer, who 
should have the means of stopping the computations and making appropriate 
changes in the problem formulation. 

Structural optimization algorithms, or more specifically, the methods 
of material redistribution, can be divided into two broad categories, which 
we call ’’direct" and "indirect” techniques. In the first category, the 
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redesign is viewed as a purely mathematical programming problem, where the 
merit function (weight) is to be minimized within certain constraints. Typi- 
cally, the heart of the method is a systematic search procedure that works 
directly with the merit function and the constraints and converges to either 
the global or local minimum weight design. 

The direct methods trace their origin to the operations research and 
optimal control theories and have dominated the literature on automated 
design since the birth of the subject, just over a decade ago. Variations of 
the technique are numerous — names like "feasible directions," "steepest 
descent, " and "gradient search" are just a small sampling of the terminol- 
ogy used. A very readable review of direct methods has been presented by 
Schmit [ 1] . 

Despite the successful treatment of many problems, direct optimization 
does not appear to be a practical method for the design of large structures. 

It has become apparent that the number of design iterations required for con- 
vergence to optimal design increases very rapidly with the number of elements. 
Present estimates [ 2] place an upper limit of 150 design variables (structural 
elements) on the size of the structure that can be treated with acceptable 
computational economy. For this reason alone, direct methods have been 
increasingly overlooked in recent work on practical aerospace structures and 
will not be discussed further in this report. 

A somewhat different approach to optimal structural design was devel- 
oped in the late 1960's. The idea was to avoid the inconvenience and compu- 
tational inefficiency of working directly with the merit function and the 
original constraints — hence the name, "indirect" method. 

Indirect design methods are, loosely speaking, counterparts of the 
variational methods used in analysis. Each method is built around an opti- 
mality criterion that serves the same function as any well-known variational 
principle of mechanics; see Reference 3 for simple examples on the derivation 
and use of optimality criteria. The analogy is accentuated by the fact that 
the optimality criteria frequently specify the energy distribution in an optimal 
weight structure. 

An optimality criterion is mathematically equivalent to the design ob- 
jective and the constraints; consequently, its use also leads to a true (local) 
minimum weight structure. The advantage of the indirect method stems from 
the presence of the behavioral variables, rather than the total structural weight, 
in the optimality criterion. In contrast, the direct method displayed the 
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behavioral variables in the constraint conditions only. Since the constraints 
are the source of most mathematical difficulties in optimization problems, the 
use of the optimality criterion usually leads to a more efficient design algorithm. 

Despite the aformentioned superiority, a rigorous application of the 
optimality criteria to realistic design problems may still require unacceptably 
long computer runs. Consequently, it is a common practice to introduce 
convenient ad hoc modifications into the formulation of the problem, which 
have the effect of sacrificing some of the structural efficiency for computational 
economy. Although the resulting design will not have the true optimal weight, 
a modified method is still a substantial improvement over the traditional design 
technique. 

Most of this report is devoted to the derivation of optimality criteria 
and their utilization in design algorithms. Although much of the material is 
derived from published descriptions of existing optimization programs, the 
report should not be considered merely as a literature survey. The publica- 
tions listed either treat specific aspects of optimization, introduce special 
techniques, or are applicable to only a restricted class of elements. In 
contrast, the present work attempts to establish a general viewpoint to struc- 
tural optimization. The optimality criteria and redesign techniques are first 
introduced in general terms, then specialized to different constraint conditions 
and compared with the methods proposed in the references. 

It is important to keep in mind that computer-oriented structural opti- 
mization is still in the developmental stage. The only existing program with a 
sufficient capacity (3000 elements, 6000 degrees of freedom) to design a large 
aerospace structure [4] is largely limited to stress constraints 2 and restricted 
to elements with special properties. The remaining programs listed in the 
references are small demonstration algorithms involving no more than a few 
hundred elements with simple size-stiffness relations. Therefore, the design 
methods proposed here and in previous publications represent only the first 
steps toward a practical, fully automated structural optimization algorithm. 


BASIC CONCEPTS 


Restrictions 

This report is restricted to the minimum weight design of linear, 
elastic structures. The layout of members and the loading (static, unless 


2. Displacement constraints are treated in an indirect, rather inefficient 
manner — see the section, A Selected Survey of Optimization Programs. 
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specified otherwise) are assumed to be given, the only design variables being 
the sizes of structural members — the cross-sectional ared and thickness for 
one- and two-dimensional elements, respectively. The size is presumed to be 
constant within each individual element. It is also assumed that finite element 
displacement methods are used in the analysis cycle. 


Constraints 

The behavioral constraints treated in the report are upper limits on 
stresses and displacements, and lower bounds on general buckling loads, 
natural frequencies, and flutter speeds. As these constraints play a major 
role in the derivation of optimality criteria, they can be classified as primary 
constraints . 

In most design problems it is desirable to incorporate additional con- 
ditions, called secondary constraints , in the design algorithm. Minimum and 
maximum limits on member sizes, prohibition of local buckling, and equal 
size constraints (the requirement that the sizes of certain members be the 
same), fall into the last category. The secondary constraints mainly deter- 
mine the details of programming, having little effect on the optimality criterion 
itself. 


An insight into the relationship between the constraints and the optimal 
design can be obtained only through geometrical abstraction. To this end, we 
introduce the concept of design space — an N-dimensional Euclidean space, 
where N is the number of independent design variables. The coordinates are 
the design variables A^ , i = 1, 2, . . . N , so that each point in the space 

represents a specific design of the structure. 

The points that do not violate any constraints are known as feasible 
designs . The boundary between the feasible region and the remaining space 
is called the surface of active constraints , and the points on that surface are 
termed as critical designs . The optimal design problem involves finding that 
point of the feasible region that is associated with the lowest structural weight; 
it is invariably a critical design. 

The three-bar truss in Figure la lends itself to a simple example of 
design space. All the members are assumed to be made of aluminum, with 
E = 6. 8948 x 10 10 N/m 2 (10 7 psi) (E is Young’s modulus). The loading con- 
sists of two alternate forces of 88. 96 N (20 kips) each, as shown. The 
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structural volume (weight) is to be minimized subject to the following con- 
straints: 

|<r | * 137. 88x 10 5 N/m 2 (20 ksi) , 

u — 0. 381 cm (0. 15 in. ) , 

v 

u^ — 0. 762 cm (0. 3 in. ) , 

and 

A. ^ 0. 508 cm (0. 2 sq. in. ) , 

th 

where cr is the stress in the i member, A is the cross-sectional area, 
i i 

and u, and u represent the displacement components defined in Figure lb. 
n v 


P = 88.96 N 


/ = 2.54 m 


I— / — )— / -| 



a. 


b. 


Figure 1. Three-bar truss. 
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The symmetry of loads and structural layout enable us to simplify the 
problem somewhat by considering only one of the loads and imposing the addi- 
tional condition A 3 = Aj (Fig. lb). This reduces the dimensions of the 
design space from three to two. 

In this simple example it is easy to analyze the structure with arbitrary 
values of the design variables. The results are 


A 2 + nT^A! 

P_ 

2A 2 + nT^Aj 

A l 

\T?a 2 

P_ 

2A 2 + \f2'A 1 

A 2 

A 2 

P 

2A 2 + nTS!’ Aj 

Ai 

Pi 



U v (\T2! A 2 + Aj)E 


9 


9 


and 




Pi 

A,E 


All the constraints can now be shown in the design space (Aj - A 2 plane) by 
plotting the lines = 137. 88 x 10 5 N/m 2 (20 ksi) , u y = 0.381 cm (0.15 in. ), 

etc. , as has been done in Figure 2. The active constraints are determined by 
the envelope of all the constraint lines. 

Figure 2 also shows constant volume contours of the structure, obtain- 
able from V = (2\T~2 A i + A 2 )i (V denotes the material volume). The 
optimal design can readily be found by inspection; it is represented by the 
point where a volume contour is tangent to the active constraint line. 
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(cm 2 ) 


CONSTRAINTS 
ACTIVE CONSTRAINTS 
VOLUME CONTOURS 



Figure 2. Design space for three-bar truss 



The location of the optimal point in Figure 2 is atypical. In the major- 
ity of problems it is found at the intersection of two or more constraint surfaces, 
as shown in Figure 3a; that is, the optimal design is commonly determined by 
several constraints simultaneously, rather than by a stationary point (Fig. 3b). 



a. Intersection of b. Stationary point, 

constraint surfaces. 


Figure 3. Examples of optimal points. 

An important feature of optimal design topology, which is apparent in 
Figure 2, is that only some of the constraints imposed on the problem are 
active. Moreover, only a portion of these constraints are active at the optimal 
design point itself, i. e. , are directly involved in determining the minimum 
weight point. 

It turns out that it is not too difficult to construct an optimal design 
algorithm if the constraints that are active at the optimal point are known 
beforehand. Unfortunately, this can be done only in a few, small-scale prob- 
lems (see Reference 5 for an example). In the majority of structures a large 
portion of the computational effort must be expended, directly or indirectly, on 
finding the constraints associated with the optimal point. 

Optimality Criteria 

Consider at first one load condition and a single primary constraint, 

q - q* , (!) 

r r 
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where q^ is a behavioral variable and q* its prescribed upper limit. In 
addition, we permit constraints on element sizes: 

.) . — A. ^ (A.) 

1 mm i 1 max 



Let A represent a design which is assumed to be critical, i. e. , 



( 3 ) 


but not necessarily an optimal one. We presume that it is possible to compute 
the derivatives 



( 4 ) 


The changes Aq of the constraint variable due to changes AA of the design 
r i 

can be estimated from 


Aq = Y. Q . AA. . (5) 

r “ n i 
i 


Equation (5) is an approximation for finite values of AA. ; it is exact only when 
|AA| is infinitesimal. 

I ^ I 

The corresponding change of the structural weight is given by 


AW = 



( 6 ) 
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where p^ is the unit weigh of the i element and represents the length 
or surface area for a one- or two-dimensional element, respectively. 

The design A can be improved if AA^ exist that do not violate the size 

constraints and for which AW < 0 , Aq — 0 . We will now investigate the 

r 

possibilities of weight reduction in detail. 

1. If s o for an element, it is evident by inspection of equations 

(5) and (6) that we can obtain a lighter-weight feasible design by simply de- 
creasing the size of that element. It follows that at optimum weight design we 
must have 


A. = (A.) . if Q . > 0 

1 1 min n 


(?) 


2. Consider any two elements, denoted by i and j , for which 
< 0 and Q^. < 0 . It might be possible to save weight by removing 

some material from element i and adding a portion of that material to element 
j, or vice versa. To avoid violating the primary constraint, we consider only 
design changes which lead to Aq = 0 . According to equation (5), the last 
constraint is equivalent to 



Q . 

ri . . 
7— AA. 
Q . 1 

rj 


( 8 ) 


The weight is reduced if AW < 0 , i. e. , 


p L AA + pLAA. <0 
i i i J J J 


Substituting from equation (8), we get 


11 


( 9 ) 




< 0 


If the size constraints are not active, we can satisfy equation (9) by choosing 
the appropriate sign for AA. as long as the term in the parenthesis is not 

zero. Consequently, the design can be an optimal one only when the term 
vanishes, i. e. , 


P.L. 
3 3 




The last expression is equivalent to the equations, 


P.L + A. Q . — 0 (10a) 

l x r ri 


and 


P.L + \ Q — 0 , (lOb) 

3 3 r rj 


where A. is a positive constant (Lagrangian multiplier); its value is deter- 
mined from the condition q = q* . The positive sign requirement follows 

r r 

the inequalities p.L >0 and Q . < 0 . 

i i ri 

Next consider the case A = (A ) . . Because of the minimum size 

i i min 

constraint, the design change is restricted to AA^ > 0. It follows from an 

inspection of equation (9) that the weight can be lowered only if the quantity in 

No weight reduction is possible if and only if 


0 . 


the parenthesis is negative. 


Q . 

p.L. - p.L. S 

3 ] Q • 
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Assuming A. to be already at its optimal value, we can use equation (10b) to 
substitute p. L./Q r - = " » obtaining the optimality criterion 


pL + \ Q . > 0 if A 
i i r ri 1 


(A.) . 

1 min 


Similarly, it can be shown that 

p L + X Q . - 0 if A = (A ) 
l i r ri i i max 


(lla) 


(lib) 


at the optimal design. 

The optimality criteria, equations (10) and (ll), can be generalized 
for cases where the optimal point is determined by several primary constraints 
and load conditions [ 6]. The results can be stated as follows: The design A 
is an optimal one if 


p L + 
l i 



0 

if 

(A.) . < A. < 

l min i 

< A i> 

0 

if 

A. = (A.) . 
i i mm 


0 

if 

A. = (A.) 
l l max 



( 12 ) 


where 


X 

r 


= 0 if q < q* 
r r 


>0 if q «= q* 
r r 


(13) 


Equation (13) requires additional explanation. If there are M load 

L 

conditions, each with constraints, the total number of primary constraints 
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in the problem is M = M T M ; consequently r = 1, 2, ... M. Only some 

L c 

of these constraints are active (q = q* ) at optimal design and thus enter 
the optimality criterion. The inactive constraints (q^ < q*) are eliminated 
by setting \ = 0 for the appropriate values of r . 

Equation (5), from which the optimality criteria were derived, is a 
linear approximation of the constraint surface in the neighborhood of the design 
A . It follows that the optimality criteria are valid only in a small region 
around A ; that is, they are conditions for local optimality and do not guaran- 
tee that A is a global minimum weight design. Proof of global optimality is 
a difficult problem, which has been resolved only for some special cases. 


Weight Reduction Cycle 

We assume again, for the sake of simplicity, that the optimal design is 

( J/) 

governed by only one primary constraint : ^ — q* . Let A be a point in 

the design space, which we call the current design. It does not have to be a 
critical design, nor does it have to lie in the feasible region. Our task is to 

modify the current design in such a way that the new design A^ + ^ is closer 

( p) 

to the optimal point than A . In particular, we want a repeated application 
of the weight reduction cycle to produce a sequence of designs that converges 
uniformly to the minimum weight configuration, as shown in Figure 4. 


A computationally efficient redesign equation can be obtained from the 
optimality criteria (12) and (13). For a single constraint, the optimality 
criterion for active members , that is, members governed by primary con- 
straints rather than size limits, is 


p i L i 


+ \ Q . 
r ri 


0 . 


Multiplying both sides of the equation by (l - a)A , where a is a constant 

i 

to be determined later, and rearranging the terms, we can write 

A. = C.A. (14) 

ill 
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Figure 4. Example of a convergent sequence of designs produced by 
repeated application of weight reduction operation. 


where 


C 

i 


a 


{1 - a) \ r 



(15) 


Equation (14), being the optimality criterion, is satisfied identically 
for each active member when A is the optimal design, regardless of the value 
of a . On the other hand, if A is a nonoptimal point, equation (14) can be 


used as the redesign formula for active members: A 


MD 


= C. A. 

l l 




Repeated applications of the formula are equivalent to the solution of the 
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optimality criteria by the method of relaxation. The parameter a , - 1 < a < 1, 
is known as the relaxation factor ; its value determines the rate of convergence. 
It has been found that the optimal rate of convergence usually requires some 
underrelaxation (a > 0). 

In order to include passive members — members governed by size 
constraints — into the redesign equation, the latter is rewritten as 


(v + 1) _ 


C.A M 
1 1 

if 

(A.) < C.A. < 

i mm i i 

^Vmax 

(A.) . 
v l min 

if 

C.A.^ < (A.) . 

li i min 

(16) 

(A.) 

v l max 

if 

C.A. ^ > (A.) 
ii i max 



where C is given by equation (15). 
i 


Before C. can be evaluated, however, the value of the Lagrangian 
multiplier must be found. As pointed out before, is determined in 

such a way that the design + is critical, i. e. , from the condition 

{v + l) 

q = q* • 

r r 


We assume for the time being that the identity of the active and passive 
members is known a priori. The change in q^ due to the design changes in 

the passive members is, according to equation (5), 


(Aq ) = T. Q . 

r pass . u „ n 
i pass 1 


(A.) . - A. 

l min l 


(4 


l Q. 


i pass 2 


(A.) - A. 

i max i 


00 ' 


(17) 
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The two sums represent the effects of members governed by minimum and 
maximum size constraints, respectively. 


The contribution of active members to change in the constraint variable 


(Aq r ) aot - 2 Q rl AA| . 

i act 


From equation (16) we obtain for active members 


AA = A.^ + 1 ^ _ A ^ = (C. -1)A.^ 
i 1 i ii 


Substituting for C . from equation (15) yields 


AA = - (1-aO (l + X %-)A M 
i l rp.1. i 


Therefore, equation (18) becomes 


(Aq ) = - (1 - a) Z VV W 

r act , u \ n r p.L. / i 

i act \ ii/ 


The total change in q is 

r 


Aq r = <Aq r>pass + (Aq Act 


We set Aq^ = q* - q^' , which is equivalent to the requirement 
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q + = q* t and substitute for (Aq ) from equation (19). The 

r r r 3-CC 

resulting equation can be solved for the Lagrangian multiplier: 


X 

r 


q - q ^ - (Aq ) + (1 - O') Yj q . a. 

r r r pass . LJ , n 1 


(p) 


i act 


Q 2 , 


<*-“> A FIT a, 


M 


i act i i 


( 21 ) 


Equation (21 ) can be evaluated only if the identities of active and 
passive members are known beforehand. Since this is generally not the case, 
the weight minimization cycle itself must be carried out by a trial-and-error 
procedure outlined as follows. 


1. Analyze the current design 
of the constraint variable. 



and compute the gradients 



2. Set A. 


(v + 1) 


bers to be passive. 


(A ) . if Q . > 0 and consider these mem- 

i min n 


3. For the remaining members, use the same division into passive 
and active groups that occurred at the end of the previous redesign cycle. 

If the present redesign operation is the first one, assume that all these mem- 
bers are active. 


4. Compute X from equation (21). 


5. 

results to 


Use equations (15) and (16) to calculate A. + ^ and use the 
reclassify the members into passive and active groups. 


6. If the classification of the members has remained unchanged, the 
redesign cycle is completed; otherwise, use the new classification to repeat 
steps 4, 5, and 6. 


Experience reported in References 2 and 7 and in a manuscript now 
being prepared for publication 3 has shown the method to be efficient. In most 

3. J. Kiusalaas, Optimal Design of Structures with Buckling Constraints. 
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cases only a few iterations are needed for the active-passive classification of 
members. Moreover, the optimization technique itself — repeated application 
of the redesign equation (16) — appears to be the most economical method of 
minimum weight design in use at the present time. Three or four cycles 
usually yield a structural weight that differs only by a few percent from the 
true optimal weight, regardless of the size of the structure. 

In the preceding discussion it was assumed that the design is deter- 
mined by a single primary constraint. It is not difficult to revise the redesign 
equations for multiple constraints (or several load conditions). Equation (16) 
will remain valid but now 


C = a 
i 


1 - a 

»ih 


E 


p act 


\ Q . 
P Pi 


( 22 ) 


where the sum is taken over all the constraints that are active at the optimal 
design. The requirement that q^ V + ^ = q* , where r applies to active 
constraints only, yields a set of simultaneous equations for \ : 


- (l - a) 


E 

i act 


Q . A 
ri i 

p i L i 


<-) 


E 

p act 



= (1 - a) ?. Q . A.^ + q* - q 
' . L ri i r j 

i act 


O') 


(Aq ) r : 
r pass 


1, 2, .. . M . 
(23) 


The contribution of passive members 
equation (17). 


(Aq ) 
r pass 


is obtained again from 


As can be seen, optimization under multiple constraints is considerably 
more complex than design with respect to a single constraint. The main diffi- 
culty is that the active constraints can be identified only by a trial-and-error 
routine, similar to the one used in separating active and passive elements. 
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The general idea is to initially assume that all constraints are active and to 
solve equation (23) for the Lagrangian multipliers. The constraints associated 
with negative values of X^ are then considered passive [ the multipliers must 

not be negative; see equation (13)] and the corresponding rows and columns 
are eliminated from equation (23). The equations are solved again and the 
process is repeated until all the remaining values of X are positive. Since the 

identity of passive and active members is also unknown at the beginning, the 
calculation for the Lagrangian multipliers must be combined with the procedure 
for classifying the elements. 

If only a few constraints are imposed, the method is practical and effi- 
cient. It has been successfully applied to design with respect to buckling, 4 
where two primary constraints (lower bounds on the first two buckling loads) 
were employed. Most practical problems, unfortunately, involve numerous 
constraints and load conditions, so that the multiconstraint optimization method 
just outlined would be exceedingly expensive. Within the confines of the present 
state of the art, economically acceptable computer times can be obtained only 
by resorting to approximate methods. Two of these methods are described in 
this report in the section, A Selected Survey of Optimization Programs. 


Behavior Modification Cycle 

It is sometimes necessary to modify the behavior of the structure by 
means other than the weight minimization operation. The sole function of the 
behavior modification cycle is to bring the design as close to the active con- 
straint surface as possible within the limitations of the linear approximation 
(5). Therefore, the criterion in choosing the appropriate redesign formula is 
the accuracy with which the behavioral changes can be predicted; the resulting 
change in structural weight is of secondary importance. 

The need for a behavior modification cycle stems from the approximate 
nature of the expression for changes in the behavioral variables (5). Exper- 
ience has shown that equation (5) is sometimes a very poor approximation if 
the design change AA is obtained from the weight minimization operation 
described in the last section. In certain cases, such as structures on elastic 
supports, 5 a repeated application of the weight minimization equations may even 
lead to divergence of th design from the active constraint surface. 


4. Ibid. 

5. Ibid. 
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In a behavior modification operation the member sizes are changed 
more or less uniformly, in which case equation (5) has been observed to be 
a more satisfactory approximation to behavioral changes. 

Figure 5 illustrates the application of behavior modification. A band 
of acceptable, near-critical designs is defined by 


e < q < q* + 
1 r r 


e 


2 


where and € 2 are predetermined constraints. Whenever the current 
design is outside the acceptable band, a behavior modification cycle is used to 
bring the design closer to the active constraint surface. The weight minimiza- 
tion formula is used only on designs that lie inside the acceptable band. 



Figure 5. Example of a sequence of acceptable designs produced by 
behavior modification and weight reduction operations. 
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This method of design, has two advantages over the repeated application 
of the weight minimization cycle alone. The first and most important of these 
is the ’’stabilizing" effect that the behavior modification has on convergence of 
the design procedure. The use of behavior modification can not only turn a 
divergent process into a convergent one but will also improve convergence in 
general. Second, the sequence of acceptable designs (points 1, 2, and 3 in 
Fig. 5) is useful in monitoring the design process. A weight comparison of 
successive designs can be used to terminate the design process whenever the 
weight reduction becomes small or ceases altogether. This method of termi- 
nation is particularly valuable when approximate, multiconstraint optimization 
techniques are used which converge to a point other than the true minimum 
weight design. 

Uniform scaling is perhaps the most efficient means of obtaining a near- 
critical design since it yields a good prediction of behavioral changes. The 
sizes of all active members are scaled by the same factor , i. e. , by the 
operation 


\v + 1) 



(24) 


The passive members can be accounted for by using equation (16) as the scal- 
ing equations, with 


C. = 

l 




(25a) 


The value of p is 
(v 

requirement q^ 


computed in the same manner as 

+ ^ = q* . The result is 
r 


\ 

r 


namely from the 




q* - q - (Aq ) + £ Q • A - 

r r r pass . u L n i 

i act 


Tj Q • A. 

. u L ri l 

i act 


(*0 


(25b) 
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In all other respects the sealing operation is identical to the iterative weight 
reduction cycle described previously in the subsection. Weight Reduction 
Cycle. 


In multiconstraint design problems a value of p should be evaluated 
for each constraint and the largest value used as the scale factor. In practice 
it is sufficient to confine the calculations to a few constraints, the choice being 
based on the magnitude of the ratio 


R = 
r 



(26) 


The constraints yielding the largest values of R are more likely to deter- 
mine the critical design. r 

Another behavior modification operation that has been used is the 
gradient travel mode [8, 9, 10] . The idea is to bring the design to the active 
constraint surface with a smaller change in weight than in uniform scaling. 

In gradient travel the weight change of each active member is proportional to 
its effectiveness in changing the constraint variable: 


9 <1 




(27) 


where p is the constant of proportionality. Substituting 


AW. = p.L (a ^ + _ A ^ 

i i i \ i i 


and 


9q /aw. = Q ./(p.L.) , 

r i n i i 


we can rewrite equation (27) in the form 
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where 


A 


(p + 1) 


= C. A 




c. 

i 


Q. 


l + ij. 


ri 


(P.L. ) 2 A.* 1 '* 
.11 1 


(28a) 


As in the case of uniform scaling, we can include passive members in 
the gradient travel mode by adopting equations (16) as the redesign formulas. 
The scale factors C. would, of course, be computed from equation (28a). 


The requirement 


O' + l) 


q* yields the value of p : 


M 



(Aq ) 

r pass 


Z 

i act 



(28b) 


The gradient travel mode represents, within the linear approximation, 

( v) 

the shortest distance in the design space between the point W and the active 

portion of the q = q* surface. In the absence of passive members the rede- 
T T J / \ 

sign vector is normal to the q = constant surface at W , as shown in 

r ~ 

Figure 6b. In comparison, the redesign vector used in uniform scaling passes 
through the origin of the design space (Fig. 6a). 


The gradient travel operation does make sense only if the behavior mod- 
ification cycle requires an increase in stiffness (and in structural weight). If 
a decrease in stiffness is required, the gradient travel mode conflicts with the 
design objective because it reduces the stiffness with a minimal saving of 
structural weight. 
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Wi = P 1 A, 


W 1 "Pl^l A i 


a. Uniform scaling. 


b. Gradient travel. 



c. Weight gradient travel. 


Figure 6. Examples of behavior modification operations. 






A behavior modification operation that is better suited for reducing 
structural stiffness is the weight gradient travel operation [ 1 0] , 


A A 9W 

AA. = n — — 
i 3 A. 


= MP L. 
l l 


(29) 


where i applies to the active members only, 
equation (29) becomes 


Using AA = 
i 


(v + 1) . M 


-A. 


A (-» 


C A. 
1 1 




where 


C 

i 


P.L. 
i i 

P — — - 1 


O') 


(30a) 


The expression for (i , obtained again from q^ + ^ = q* , is 


P = 


q* - q ^ - (Aq ) 
r r v r pass 

E P.L. Q . 

, u L i i n 
i act 


(30b) 


The weight gradient travel mode in the design space is shown in Figure 6c. 
The design change vector is normal to the constant weight contour, the direc- 
tion of steepest descent on the weight surface. 
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GRADIENTS OF CONSTRAINT VARIABLES 


General Considerations 

The key to indirect optimization (as well as direct optimization) lies 
in the computation of the constraint variable gradients. The calculations are, 
with the exception of stress constraints, straightforward and economically 
feasible, provided that the relationships between the member sizes and the 
element stiffness matrices are specified in advance. 

We denote the stiffness matrix of the i^ element by [Kj and the total 

stiffness matrix of the structure by [ K] . The generalized displacement vec- 
til 

tors of the i element and the structure are written as { u.) and { u) , 

respectively. The rules for compiling [K] from the element stiffness are 
based, as usual, on the invariance of the total strain energy U : 

U = |{u} T [K]{u) = ^{u.lV.Uu.} + U g 


where the sum applies to all the elements in the structure and U denotes the 

b 

strain energy of elastic supports. 

To account for equal area constraints, we introduce the group stiffness 
matrix [ K ] and the group displacement vector { u ) , defined by 


U 

g 


} T tK ){u } - Z {“/tK.Uu,} 

6 b b j group g 


where the sum is taken over the elements belonging to the g equal size 
group; that is, elements that must have the same size A . It is assumed 

that all the elements in any one group have an identical size-stiffness relation. 
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Henceforth, no distinction will be made between problems with equal 
size constraints and those without the constraints. All the formulas that follow 

are valid for both types of problems, provided that A. , { u. } , [ K. ] , etc. , are 

th 111 

interpreted as belonging to the i group if equal size constraints are applied. 

In other words, all members belonging to an equal size group are simply lumped 
together into an equivalent element for redesign purposes. 

Development of the general theory requires no a priori knowledge of the 
exact form of the element size-stiffness relations. It is sufficient to assume 
that the derivatives of the stiffness matrices, namely 3[Kj / 9A^ , are cal- 
culable. On the other hand, the size-stiffness relations play an important role 
in some programming details, particularly in the storage of the element stiff- 
ness matrices. 

The great majority of design problems can be handled by using element 
stiffness matrices of the form 


ik,] - l 

m= 0 


k. 

1 


(m) 


A m 
A i » 


(31) 


where the unit stiffness matrices 


(m) 


are independent of the element 


—4 

size A^ . The stiffness due to the nonvariable portion of the element is repre- 
by and the direct stresses contribute to k.^ 1 ^ 


sented 


bending and twisting stiffnesses are included in 


1 


(3) 


, (1) ]a., 


A ; the 

l 


k. 

l i 


( 2 ) 


or 


A^ , depending on the type of the element. Three samples of size- 


stiffness relations are shown in Figure 7. 

The derivatives of element stiffness matrices are 


[K. .] = 
1,1 


a[K ] 

i_ 

a a. 


- 2 , » W m) ] ■ 


(32) 
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CHANGE IN GAGE THICKNESS t ONLY: 

l x - c 1 A I y — c 2 A J 2 c 3 A 3 


PROPORTIONAL CHANGE IN ALL DIMENSIONS: 
l x = c 4 A 2 ly = c g A 2 J = c 6 A 2 



a. 

CHANGE IN THICKNESS OF REINFORCEMENT t 
ONLY: 


lx~ c o + c i A ly - c 2 + c 3 A 

J = c 4 + c 5 A + c 6 A 2 + c 7 A 3 


b. 



Figure 7. Examples of size- stiffness relations. 
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Since both [K. ] and [ K. . ] are used once in each design cycle (the stiffness 
1 1,1 

matrix is required in analysis and its derivative in the redesign cycle), it is 

necessary to store £k. , m = 0, 1, 2, 3 separately for each element. 

This would require a modification of conventional analysis programs if they 
are to be used in an optimization algorithm. 


A special case is obtained when 


[K. ] = 

l 


t . (m) ] A i m 


(33) 


The size- stiffness relation in equation (33) makes the uniform scaling opera- 
tion particularly simple if all the members are active (no size constraints). 

(u + l) (v) 

If p is the uniform scale factor, such that = p A , then 


and (34) 


where {f} , (u) , and p represent generalized forces, displacements, and 
buckling loads, respectively (it is assumed that applied loads do not depend on 
member sizes). Consequently, no reanalysis is required after a uniform 
scaling cycle. 

The linear size-stiffness relationship (m = l) is of particular signifi- 
cance in aerospace structures. Shear panels, membrane elements, and thin- 
walled beams with variable gage thickness (Fig. 7a) fall into this category. 
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Uniform scaling, together with equations (34), has been extensively used in 
aerospace-oriented optimization programs [2, 4, 7, 11, 12, 13]. 

Although the simple uniform scaling procedure (scaling of all members) 
appears to be very attractive, its value is somewhat dubious because the opera- 
tion is inconsistent with constraints on member sizes. The selective scaling 
procedure in equation (16), based on active-passive member categories, seems 
to be preferable. 


Stress Constraints 


Rigorous design with stress constraints is still an unresolved problem 
in the optimization of large structures. The approximate methods that have 
been used, and are described in this report, fall somewhat outside the general 
optimization theory developed in the preceding sections. 

Because of notational complexities, it is difficult to present the ideas 
in any generality. They are best displayed for the simplest class of problems 
— the design of trusses. 

We consider first the constraints on the tensile stresses, 


< 7 — cr* 

r r tens 


The last inequality is identical to equation (l) if we set a 


q and 
r 


<T* 

r tens 


q* . For members of a truss, 
r 


a = 


r A 


where is the axial force in member r . Consequently, 


da 


\(v) 


Q - = 


ri \ 3 A 




(v) 

r - 1 

2 6 ri + 


3P 


» 


A {v) V A i 

r 


(35) 
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with 6 . denoting the Kronecker delta, 
ri 

It is a common practice to neglect the last term in equation (35), i.e. , 
to assume that the change in member forces during redesign is negligible . 

This assumption leads to an enormous simplification of the design algorithm, 
apart from avoiding the difficulties of computing 9P^ / 3 A. . 

The advantages of the simplified formula, which can be written 


Q ri a ^ 6ri 


(36) 



are most apparent in multiconstraint design problems. It can be verified by 
an inspection of equation (23) that the equations for the Lagrangian multipliers 

A become decoupled upon substitution of equation (36). As a result each 
P 

stress constraint and load condition combination can be handled as a single- 
constraint problem. The largest member size predicted by all the load condi- 


tions is chosen as the new design 


(v + 1) 


Substituting equation (36) in the equations of the weight reduction cycle, 
we obtain from equation (21), 


A 

r 


(2 - a ) cr 


(v) 


- cr 


r tens 


(1 - «) 


(^y 


p L A 
r r 


(v) 


Equation (15) then gives the redesign scale factor for the 


th 

r 


element: 


C 

r 


a 


+ 


(2 - a) a ^ 
\ r 



- cr* , 
r tens 


if 



(37a) 
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Constraints on compressive stresses 


cr — - cr* 
r r comp 


can be handled by setting q = - <r in equations (21 ) and (15). The result 
is r 


C 

r 


a + 


(2 - a) <t ^ + cr* 

r r comp 

„ (v) 


if 



(37b) 


As mentioned before, if there are several load conditions, a value of 
is computed for each condition and the largest value is chosen as the scale 

factor in the redesign equations: 


( 

C.A.^ 
1 1 1 

if 

C A > 

i i 

(A.) . 
l min 


i ' 

) 

] 




(38) 

( 

' < A i>min 

if 

C.A. (V) < 

l l 

(A.) . 
l min 

• 

The statement A 

i 

+ 1) = (A.) 

max 

if C.A.^ 

l l 

> (A.) 

l max 

has been 


omitted from equation (38) since the simplified design equations are obviously 
incapable of handling upper limits on the design variables. This is one of the 
major drawbacks of the approximate method. 


It can be seen from equations (37a) and (37b) that no change in the 

(v) _ _* nr , ^ (v) 


member size occurs 


(C = l), when a = a * . or cr ' = - a* 

r r r tens r r comp 


Therefore, the final design is known as the fully stressed design, where each 
element is stressed to the allowable value under at least one load condition or 
is governed by the minimum size constraint. 


The fully stressed design always coincides with the minimum weight 
design for statically determinate structures. In the case of static indeterminacy. 
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it is generally just an approximation to the optimal weight configuration. The 
error (weight penalty) is claimed to be small in most practical design prob- 
lems but this contention has not yet been proved conclusively. The fact is that 
special problems can easily be conceived where the discrepancy between the 
minimum weight and the fully stressed designs is significant. 

Equations (37a) and (37b) have actually not been used in the existing 
optimization programs. The fully stressed design is traditionally obtained 
more directly by the stress ratio redesign method [2, 4, 7 ] . The assumption 
that member forces do not change with redesign is equivalent to the equations, 


„ A M = a* A (l, + 1) 
r r r tens r 


if > 0 

r 


and 


- <r M A (,,) = <r» A 

r r r comp r 


( I ' + 1 ) „(") 


if a v 7 <0 

r 


Solving the equations for A ^ + ^ yields the scale factors to be used in the 
redesign formulas (38): 


C = 
r 


<J ^/cr* 
r r tens 


- a ^ /a* 

r r comp 


if 


if 


> 0 
r 


cr^ < 0 
r 


(39) 


No benefit has been obtained by introducing a relaxation factor into the stress 
ratio redesign equations [ 14 ] . 

For a general finite element the stress vector at point P can be 
obtained from the nodal forces {f } by the operation 
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( 40 ) 


R!-[<]<v 

r p] 

where the matrix S. depends on A. and the location of P . The stress 
constraint 


°i eff 



( 41 ) 


is imposed on the effective stress in the member, the latter being given by 


i eff 



( 42 ) 


The operation T is determined by the failure criterion and is independent of 

A. and the location of P . 
x 

For example, if the Von Mises yield criterion in two dimensions is 

used, 


°i eff 


= ” ax + (°i P )y - ( ff i P )x(°i P ). + 3 ( T i P ) 


P \ 2 
xy 


The stress ratio redesign is usable only for the special case 



( 43 ) 


where ^ 
to obtain 


cr 


does not depend 

(v+ 1 ) 


i eff 


= oy 


on the member size. 
, we must use 


It is readily verified that 
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( 44 ) 


C. 

1 



1 

P 


in the redesign equation (38). All of the elements used in References 2, 4, 7, 
11, and 12, apart from having linear size- stiffness relations, are of the type 
described by equation (43), with /3 = 1 . The only exception is the plate ele- 
ment used in Reference 4. 

Stress-constrained design with more general elements appears to be a 
very difficult problem and thus far no entirely satisfactory methods have been 
proposed. Certain elements with a relatively simple relationship between the 
effective stress and the nodal forces, such as plates, and circular thin- walled 
tubes with a fixed t/r ratio seem to be manageable (we assume combined 
action of direct forces, bending, and torsion). However, these elements must 
be designed by equations similar to (37a), (37b), and (38) since the stress 
ratio method would not be applicable. The alternative is to use specially 
tailored failure criteria, as has been done in Reference 4. 

Apart from the failure criterion cr — a* , an optimization algor- 

i eii i 

ithm should also provide for local buckling constraints. Again, one can make 

use of the assumption that nodal forces are unchanged during redesign. With 

the nodal forces known, the minimum element size required to prevent it from 

buckling can be computed provided, of course, that the buckling design data for 

the elements are available [4] . This size is then compared with (A.) . and 

l min 

the largest value used as the lower bound in the design with respect to stress. 
The major difficulty of the method lies in providing buckling design data for 
each element. 

Before leaving the subject of fully stressed design, another flaw of the 
method should be mentioned. The finite element, stiffness method generally 
does not yield stress fields that are continuous between elements, i. e. , they 
do not satisfy equilibrium conditions exactly. As a result the member sizes of 
the final design reflect these discontinuities and may produce a weight distribu- 
tion that has an intuitively "wrong” appearance. 

The problem is solved in Reference 4 by a stress smoothing procedure, 
called the nodal stress method . Roughly speaking, the method redistributes 
the nodal forces predicted by the finite element analysis in such a manner that 
the net force acting on each node vanishes. A bonus of the nodal stress method 
is a considerably faster convergence of successive designs to the fully stressed 
state. 
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The difficulties with abnormal weight distribution have been tackled in 
References 11 and 12 by an energy approach . The scale factor in the redesign 
equations (38) is chosen as 


C. - 
i 



9 


(45) 


where IT is the strain energy of the i^ member (the maximum energy pro- 
duced by any one load system is used if multiple load conditions exist) and 
U* is called the strain energy capacity of the element. The latter is defined 

as the energy stored in the member at failure, assuming a uniform strain field 
(constant tension in burs, constant bending in beams, etc.). The constant X 

is the same for all elements and is adjusted so as to bring the design A. + ^ 

to the active constraint surface. 


If the structure is composed of constant strain members, such as bars 
and membrane elements, the energy approach is identical to the stress ratio 
method (\ = 1 would be used). A comparison of the results produced by the 
two methods is nbt available for more complex elements. 


Displacement Constraints 

A displacement constraint can be expressed as 


u ^ u* 
r r 

where u* is the prescribed upper limit on the generalized displacement u^ . 

The above inequality has the same form as the standard behavioral constraint 
(l). Since the displacement gradients can be calculated in a straightforward 
manner, the design algorithm can be constructed directly upon the general 
theory developed in this report. 

The most economical way of evaluating displacement gradients is based 
on the dummy load method . To obtain the derivatives of the generalized dis- 
placement component u^ , we first place a unit dummy load on the structure 
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in the direction of the r-coordinate. Denoting the dummy load vector by 

, where f ^ = 6 , and the resulting displacement vector by ^u^j, 

equilibrium equations of the structure under the dummy load are 


Kl{u (r) }= {f (r) } 


(46) 


T 

Multiplying both sides of the equation by (u) , where (u) is the displace- 

ment vector due to real loads, we obtain 

(u} T [K] {u (r) }= u f . (47) 

The right-hand side of equation (47) was obtained by 


(u> T {f (r) } - E “i * lr - \ 


Differentiation of equation (47) yields 



( U ) T [ K i J i]{ u<r) } + { U>i } T lM{u (r) } 

+ (u} T [K] |u , (r) | 


(48) 


where we used the notation {u .} = 9{u}/9A. . Equation (48) can be 

> 1 1 

simplified considerably. Differentiating both sides of equation (46), we get 


[K] 


{■/I --MM 
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Similarly, differentiation of the equilibrium equations due to real loads, 

T T 

[K] {u} = {f} or {u} [K] = {f} , yields, assuming {f} to be indepen- 

dent of A , 
x 


{u } T [Kl = -{u} T [K .] . 

y 1 y 1 

Consequently, equation (48) becomes 


9u 

r 

9A. 


- { u} T [K ,](u (r) } 

y 1 


Only the stiffness matrix of the i element contributes to [ K ] , i. e. 
[K .] = [K. .]. 

,i i,i 

Therefore, 




- - {<*,} 


[K. 


1.1 



(49) 


References 2 and 7, in which the dummy load approach was first used, 
imply that equation (49) is an approximation, valid only when the changes in 
internal forces are negligible during redesign. The equation is, in fact, an 
exact expression for the displacement gradients within the framework of the 
finite element theory. 


( £ ) 

The displacements {u} and {u v '} can be calculated simultaneously 

(r) 

during the analysis of the structure by adding {f' / } to the matrix of real 
loads. The extra cost of computation would be relatively small. 
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The special size-stiffness relationship (33) produces an interesting 
result. Because 


[K.l 

1 



9 


we have 


[K. .] = m [K.]/A. 
1,1 ii 


Therefore, equation (49) takes the form 



- m U. (r) /A. 
i ' i 


> 


where 



may be called the ’’dummy energy’’ of the i element — the work done by the 
real internal forces as they undergo the displacements of the dummy load. The 
optimality criterion (12) now becomes for active members 


m 



/ W 

i 


= 1 


9 


where 


W. = p. L. A 

i ill 
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element. For a single constraint, the criterion is 


is the weight of the 


.th 

1 


U. 

1 



= c , 


i. e. , the dummy energy density must be constant throughout the structure. 
The value of c is determined by the allowable displacement u* . 


Buckling Constraints 

We assume that all the loads acting on the structure can be considered 
to be proportional to a single parameter p . The values of p at buckling are 
denoted by p^ and are presumed to be arranged in an ascending order: 


Pi - P2 - P3 • • • 


If the structure is to be safe against buckling, the constraint conditions 
are 



p* 


9 


where p* is the desired value of p at buckling. The inequality can be brought 

to the standard form of equation ( 1 ) by multiplying both sides by minus one and 

setting q = - p , q* = - p* . 
r r r 

At casual glance it may appear that the design could be based only on 
p t , the fundamental buckling load. This approach is indeed adequate if 
Pi < p 2 at the optimal design. It can be shown, 6 however, that the minimum 
weight structure may possess two fundamental buckling modes (pj = p 2 ), 
which requires the use of a multiconstraint design approach. (There is a 
slight possibility of having more than two fundamental modes at the optimal 
design. ) 


6. Ibid. 
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The buckling problem is governed by the incremental equilibrium 

equation 


[ K] {u} = p [H] {u} 


( 50 ) 


where {u} is the vector of generalized buckling displacements, and [ H ] 
represents the geometric stiffness matrix of the structure. The latter is 
symmetric and is assumed to be independent of the member sizes. 7 

Differentiating both sides of equation (50) with respect to A. , we 

obtain 


[K .] {u} + [K] {u .} = [H] (u> + p[H] {u .} . (51) 

T T 

Multiplying equation (50) from the left by (u .} and equation (51) by {u} 

* 1 

and then subtracting (50) from (51) yields 


{ u} T [K ] {u} = 

9 1 


(u} T [H]{u} . 

i 


Finally, upon substituting [K .] 

9 1 


[K. .] , we obtain 




{u i ( ->} I [ K ial {u i (r) } 
|u (r) } T [H] { U (r) } 


(52) 


7. The assumption is strictly valid only when the prebuckling state is statically 
determinate. For the case of static indeterminacy, the assumption is a 
convenient approximation; optimal design can be obtained by recomputing 
the forces in the prebuckling state after each redesign cycle. 
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til 

The superscript (r) signifies that the displacements of the r buckling 
mode are to be used in the equation. There is a striking similarity between 
equations (52) and (49), which is accentuated when the special size- stiffness 
relations (33) are used. With the substitutions 


U 

i 


(r) 


(u 1 (r) } T fK 1 ,{ Ui (r) } 


» 


G 


(r) 


{ u (r) ) T [H]{ u (r) } 


9 


the optimality criterion becomes 



= 1 


( r ) 

for the multiconstraint case and U. VW. = c for a single constraint. As 

l i 

in the design for displacement constraints, the last equation also requires a 
uniform energy density — in the present case the strain energy of buckling — 
throughout the structure. 


Natural Frequency Constraints 

Constraints on natural frequencies are handled in essentially the same 
manner as buckling constraints. We introduce p^ = u ^ 2 , Pi — p 2 — P 3 • • • » 

where w^/ (27r) are the natural frequencies of the structure. It is assumed 

that the design objective is to eliminate all frequencies below a certain value 
o>* . Consequently, the behavioral constraints are, as they were for buckling, 


P 


r 


> 


p* 


Again, the standard form of the behavioral inequality is obtained with 

q = - p and q* = - p* . 
r r r 
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The importance of optimizing the design with respect to at least two 
modes, which was discussed in the preceding subsection, is also applicable 
to frequency-constrained design. 

The free vibration equation is an eigenvalue problem of the same form 
as equation (50): 


[K] {u} = p [M] {u} , (53) 


where {u} represents the buckling mode and [M] is the mass matrix. If 
the rotary inertia is neglected, the element mass matrixes can be written as 


A , (54) 


A . In addition to the contribution of the 
1 

individual elements (54), [M] is also allowed to contain nonstructural 
inertia terms (due to masses attached to the structure). 

The gradients of p^ are obtained in a manner identical to the method 

used for buckling constraints, and the derivations are not repeated here (one 
must not forget, however, that the derivatives of [M] are nonzero). The 
result is 


t M. ] = 

(o) 

m. 

+ 

l 

i 



m. 


where 


(k)l 

v ’ are independent of 
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For the special case [K.l = and [Mj = 

we introduce the potential and kinetic energies of the i^ element, 

u i (r) - | u i (r) } T[K i , { u i <r> } 


and 

T , (r> - P r {“ i (r) } T 'M i lf i (r) } ■ 

respectively, and the kinetic energy of the entire structure 


,(r) 


P r (u (r) } T [M]{u (r) } 


Equation (55) then is 


Q ri = -(p/ T<r) )( mU i <r) 




(v) 


and the optimality criterion (12) becomes for the active members, 



If a single constraint is used, this reduces to 
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Before leaving the topic, we should point out a peculiarity of the uniform 
scaling operation. Assuming the special size-stiffness relation just described 
and an absence of nonstructural inertia terms in [ M ] , the uniform scaling 
operation, when applied to all members, would result in 


[K] 


(u + l) 



9 


[M] 


(*+ 1) 



9 


where p is the scaling factor. The vibration equation (53) for the new design 
would be 


[K] 


(v){ u) (^i) 



+ 1 > 


A comparison with the equation of the previous design, 


[Kl W {u 


= p (, ' ) [Ml W {u} W 


leads us to the conclusion that 


and 


u 


(^i) _ » 


= u 


> + l) _ 


m-1 ( v ) 
M P 
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We note that when m = 1 , no change in the natural frequencies will 
occur. Consequently, the uniform scaling operation will be ineffective in modi- 
fying the behavior of the structure if the size-stiffness relations are linear 
unless the terms that are independent of dominate the mass matrix. 


Flutter Velocity Constraints 

Design with respect to flutter velocity completes the trio of eigenvalue- 
constrained problems, the other two being designs with respect to buckling and 
natural frequency constraints. Denoting the flutter speeds of the structure by 
V , Vj ^ — V 3 . . . , and the desired lowest flutter velocity by V*, the 

constraint conditions are 


V > V* 
r 

Therefore, q^ = - and q£ = - V* . 

The equations governing steady- state motions of the structure can be 
written in the form [ 10] , 


[K ] { u} = p ([Ml + [ A ] ) { u} (56) 

where [A] is the air force matrix and p = cu 2 , being the frequency 

of flutter oscillations. The air force matrix is complex, asymmetric, and a 
function of the reduced frequency 


k = 


bco 

V 


(57) 


b being the semichord length. The exact form of [A] depends on the aero- 
dynamic theory used. 
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Multiplying both sides of equation ( 57) by V and differentiating, we 


obtain 


k 


9 V 
9 A 




b 9 p bco 9 k 

2 9 A. k 9 A. 

l l 


from which 





2co k 9 A. 
r r i 


+ 


bw 9k 

r r 

, 2 9 A. 

k i 

r 


(58) 


To obtain 9p/9A. , we differentiate both sides of (56): 


JX.j{u} + [K]|u,.j = |^-([M] +[A]){u} + p([M,.] + [A,.]^{u) 

+ p([M] + [A]^|u,. | 


(59) 


The next step is to eliminate the derivative of the eigenvector { u . } . Since 

[A] is asymmetric, we need the help of the associated eigenvector {v} , 
given by the solution of 


{v} [K] = p {v} T ([M3 + [ A 


(60) 


It can be shown that the eigenvalues p^ of equations (56) and ( 60) are identical. 


T 

Multiplying equation ( 59 ) from the left by { v } and equation (60 ) 

from the right by { u } and subtracting, we have 

,i 
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{v} T [K ,I{u} = f*-{v} T ([Ml + [Al){u} 

’ i 

+ p(v} T ([M .1 + [A j){u} 

Substituting 


[A J = (9 [A] / 3k) (ak/3A. ) 


and using the notation 


3[A] /3k = [A’] 


we get 


l£_ 

9A. 


{v} T [K .]{u) - p( v} T ([ J + |~[A’]){u} 
{v} T ([M] + [A]){u] 


' (61) 


At this point we recall that if an arbitrary value of V is used in equa- 
tion (56), the resulting eigenvalues are generally complex. A real eigenvalue 
p^ , signifying a steady- state motion, can be obtained only when V = , 

r = 1, 2, 3... . Since we are designing with respect to flutter, i. e. , steady- 
state oscillations, we must use p = p^ and k = k^ (both real) in equa- 
tion (61) and also restrict 3p /9A and 9k /9A. to real values. The 
' r i r i 

last requirement essentially establishes an interdependence between 


and 



T. ( 3p / 3A. ) dA. 
V r i i 


49 



dk = y. (9k / 9A. ) dA. 
r Lj r 1 1 


which assures us that p + dp and k + dk also correspond to the 

r r r r 

flutter conditions (are real) for the design A,. + dA^ . 

Following the technique developed in Reference 10, we separate the 
terms appearing in equation (61) into real and imaginary parts: 


{v (r) } T ([K .] - P r lM ;1 l){u W } 


, , , (r) 


P r {v< r >} T [A']{u (r) }=R 2 (r) + il 2 (r) 


(62) 


and 


{v^} T (tM] + [A]){u} = + il 3 ^ 

With the terms introduced, 9p / 9A. can also be divided into real and imagi- 

r i 

nary portions: 



9A. 

l 



(R S Cr >y ♦ (i 3 <r) ) 2 


+ 



- i. 


t \ 9k 
(r) r 


jrW r) - ( R - 





i 


( Rs (r) ) 2 ♦ (% <r) ) 2 

(63) 
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The last term in equation (63) must vanish if 3p /9A. is to be real, which 
yields r 1 


( r ) ! ( r ) _ R (r) l (r) 


3k R„ 

r 3 1 13 


3A i " R (r) I <r) - R < r > I < r > 


(64) 


From the real part of equation (63) we now obtain 




( r ) _ r ( r ) T ( r ) 

3 13 



(r^r^ - i( r M r )) ( 

V 2 3 2 3 /' 

V r M r) - R< r M r) ) 

^31 13/ 

(R< r >i< r > - B WiW) 
V 3 2 2 3 / 

(*3 <r) ) 2 * (‘3 <r) f 



(65) 


Substitution of equations (64) and (65) in equation (58) completes the expression 

for Q . . 
ri 


The main difficulty with flutter optimization appears not to be in the 
redesign but in the analysis — the solution of the flutter equation for a realistic 
air force matrix [ A ] . 


A SELECTED SURVEY OF OPTIMIZATION PROGRAMS 
Stress and Displacement Constraints 

One of the most troublesome aspects of optimal design is the treatment 
of multiple constraints other than stress constraints. The difficulties are most 
acute in displacement-constrained designs since it is not unusual to have a very 
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large number of displacement limits placed on a single problem. As noted 
before, a rigorous use of optimality criteria is out of the question and approxi- 
mate techniques must be found. 

An effective method for handling stress and displacement constraints 
has been developed by Gellatly and Berke [2,7], In essence, they consider 
the stress constraints and each displacement constraint-load combination as a 
separate, autonomous optimization problem. If there are M displacement 
constraint-load combinations, then M + 1 different values are calculated for 
each element size during a redesign cycle (one value from stress constraints 
and M values from displacement constraints ) . The largest value is selected 
as the size for each element. 

The separation of members into active and passive categories during 
each redesign cycle is accomplished by successive iterations, as shown in 
Figure 8. At first all the members are assumed to be active. During subse- 
quent iterations the active members in each of the autonomous problems are 
limited to elements that were controlled by the same problem in the previous 
iteration. The procedure is repeated until no change takes place in the active- 
passive member categories. 

Following each redesign cycle, the design is analyzed and scaled uni- 
formly to the active constraint surface. Since only elements with linear size- 
stiffness relations are used in the program, the scaling operation predicts the 
corresponding behavioral changes exactly, eliminating the need for further 
analysis. The design procedure is terminated when the structural weight 
ceases to decrease between two successive critical designs. 

The algorithm of Gellatly and Berke, like all approximate methods, 
does not converge to the true minimum weight design. The weight penalty has 
not yet been evaluated. 

The redesign equation used during the weight reduction cycle in 
References 2 and 7 differs somewhat from the formula (15) proposed in this 
report. For elements with linear size-stiffness relations, equation (15) 
becomes 


C. = a + (1 - a)\ U. v 1 /W 
l r i i 


Gellatly and Berke, on the other hand, use 
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Figure 8. Flow diagram for program in References 2 and 7. 
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( 66 ) 


C. = n/a u.^/w. 

x r 1 1 


(i) 

which can be derived from the optimality criterion U. '/W = 1 upon 
the assumption that the internal forces do not change during redesign. 


It is easily shown that the discrepancy between the two expressions for 

C is negligible for small design changes AA. (i. e. , if C. « 1 ). With 
i ii 

C. ~ 1, we have 


J 


X U (r) /W. 
r i i 


it 


1 + X U 
r i 


(r) 


/w i) 


y 


which coincides with equation (15) if we take a = 


1_ 

2 * 


Since the redesign 


formulas are strictly applicable for small design changes anyway, neither 
formula can be considered to be more ’’accurate” than the other. Equation 
(15), however, has the advantage of greater flexibility: It is applicable to 
size-stiffness relations other than linear, it allows control of the convergence 
rate through the relaxation parameter a , and it can be used in rigorous 
multiconstraint design (the simultaneous equations to be solved for X ^ would 

be linear, as opposed to nonlinear equations if the Gellatly-Berke redesign 
formula were used). 


A different approach to multiconstraint design is used by Dwyer et al. 
[4] (Fig. 9). The first phase of the design algorithm consists of repeated 
applications of the stress ratio redesign formula and the uniform scaling opera- 
tion. The displacement constraints are accounted for by computing the uniform 
scale factor from all the behavioral constraints. Consequently, each design 
cycle produces a critical (usable) design. This procedure is repeated until 
the structural weight ceases to decrease. 

The second phase is a displacement- constrained design algorithm that 
is used only when the last critical design was governed by a displacement con- 
straint. The method of redesign used in this phase may be classified as a 
gradient search procedure, consisting of alternate steps of uniform scaling and 
gradient travel mode. It is not based on the optimality criterion and is not 
always effective in reducing weight. 
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Figure 9. Flow diagram for program in Reference 4. 
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The optimization program of Venkayya [11, 12] uses essentially the 
same approach as Reference 4. The differences are minor: Venkayya replaces 
the stress ratio redesign formula with the energy ratio formula (45) and uses 
a slightly modified gradient travel mode in the displacement-constrained design 
phase. The method for computing the displacement gradients is inefficient in 
Reference 4 as well as in References 11 and 12. 


Buckling Constraints 

There appear to be no published accounts of buckling- constrained 
optimization programs. The ideas developed in this report were used by the 
author to assemble a frame optimization program, 8 which works well on a 
variety of problems where the prebuckling state is statically determinate* The 
flow diagram of the program is shown in Figure 10. 

The program accepts elements with size-stiffness relations 


[ K. 1 



9 


where m = 1, 2 or 3. In addition, elastic supports (discrete nodal supports 
or uniformly distributed element supports) are permitted. Minimum permis- 
sible element sizes and equal size constraints are also included in the program. 

It was found essential to consider the problem as one of multiconstraint 
design, where the weight reduction cycle takes into account the first two buckling 
modes simultaneously [see equations (22) and (23)] . Without this feature, the 
design did not converge in cases where the first two buckling loads were equal 
at the optimal design, i. e. , if the optimal design was located at the intersection 
of pj = p* and p 2 = p* constraint surfaces, as indicated in Figure 3a. 

The redesign is carried out by either applying the weight reduction 
equations (22) and (23) or the uniform scaling operation, depending on whether 
the critical load of the current design lies in the acceptable band or not. In 
other words, the design process shown in Figure 5 is used. The program is 
stopped if the critical load is within the acceptable range and the optimality 
criterion (12) is satisfied within a prescribed latitude. 


8. Kiusalaas, op. cit. 




Figure 10. Flow diagram for program referenced in footnote 3. 
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Experience with the program indicates that the initial rate of conver- 
gence is very fast. Three or four redesigns were usually sufficient to produce 
a structural weight within a percent or two of the true optimal weight. It was 
noted that in some problems the structural behavior (buckling modes ) became 
very sensitive to small design changes near the optimal point. This caused 
successive designs to oscillate between two near-optimal points and made 
further weight reduction impossible unless the relaxation factor was increased 
close to unity. Fortunately, in all cases tested, the oscillatory behavior 
occurred only when the design weight was within a percent of the optimal weight. 


Natural Frequency Constraints 

The first procedure for optimizing a complex-structure with respect 
to frequency constraints was published by Rubin [8, 9] . The program treats 
frames and restricts elements to linear size-stiffness reductions but accepts 
nonstructural masses. 

The layout of the program, shown in Figure 11, is somewhat similar 
to the algorithm just described for buckling-constrained design. Two frequency 
modification modes are used: uniform scaling if the frequency is to be reduced 
and gradient travel [ equation (28)1 if an increase in frequency is required. 

The weight minimization mode adopted by Rubin is a numerical search 
procedure known as gradient projection search. The redesign formula is 
AA. = pg^ , where g^ is chosen so as to maximize the weight loss - AW , 

subject to constraint Ap t = 0 (no change in fundamental frequency). The 
magnitude of the redesign vector p is obtained by trial and error. 

As was already stated in the introduction, numerical search techniques 
are generally inferior to indirect design methods that are based on the opti- 
mality criteria. This is particularly true for the design of large structures. 
Another flaw of the program is the single-constraint design approach, which 
restricts its application to structures where the optimal design is a stationary 
point on the Pi = p* constraint surface (Fig. 3b). 

A more recent frequency-constrained optimization program that does 
make use of the optimality criterion has been published by Venkayya et al. 

[13] . The: elements are again restricted to linear size-stiffness relations. 

The weight minimization mode consists of the redesign equation 
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Analyze initial design 


fundamental frequency 
above the acceptable 
\ band? / 


/ Is the 

'fundamental frequency' 
\ below the acceptable / 
\ band? / 


Frequency 

modification 
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modification 

mode-gradient 
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weight decrease below 
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fundamental frequency 
within the acceptable 
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Figure 11. Flow diagram for program in References 8 and 9. 
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c A. 

i i 


A^ +1 > . 


M 

9 


where 



( 67 ) 


The origin of the formula is the optimality criterion 


(c/W.) (U. (r) 

l i 



Similarly to equation ( 66) , equation (67 ) can be shown to coincide with the 
redesign formula (15) proposed in this paper, 


C. 

i 


a + (l - a ) c 


(r) _ T (r) 


W. 

1 


9 


provided that we take n = — and consider small design changes (C. « l). 

2 i 


The algorithm in Reference 13 does not use a frequency modification 
mode — the design changes are obtained entirely by repeated application of 
equation (67). Experience with buckling- constrained design 9 has shown that 
the use of the weight minimization mode along can lead to a divergent design 
sequence in certain problems. It appears, therefore, that the absence of a 
frequency modification mode and the single- constraint approach make the pro- 
gram described in Reference 13 applicable to a limited class of problems only. 


Flutter Velocity Constraints 

Optimization with respect to flutter has been confined to very simple 
structural configurations. The most advanced paper published to date, written 


9. Ibid. 
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by Rudisill and Bhatia [ 10] , is rather similar to the frequency- constrained 
design algorithm used by Rubin [ 8, 9] and consequently suffers from the same 
drawbacks. The main contribution of the paper lies in the derivation of the 
expressions for the flutter velocity gradients. 

The weight minimization mode used by Rudisill and Bhatia is also a 
gradient projection search procedure, where the incremental increase in the 
flutter speed is maximized while the total weight is held constant. For the 
behavior modification cycle they use flutter velocity gradient travel [ equation 
(28)1 to increase the flutter velocity and weight gradient travel [ equation (30 )] 
if a decrease in the velocity is desired. 


CONCLUSIONS 


There is little doubt that computer-automated, minimum weight design 
is an eminently practical means of structural design, even in its present stage 
of development. It is safe to predict that by the end of the next decade most 
structures, in aerospace as well as civil engineering applications, will be 
computer designed. 

The principles and methods found in the present state of the art appear 
to be further advanced than their application; that is, none of the published 
design algorithms make full use of the existing knowledge and experience. Part 
of the blame must be placed on the high cost of program development — an 
optimization program requires about twice the programming effort of a corres- 
ponding analysis algorithm. In addition, structural optimization is still a 
peripheral area of structural design, known only to a small group of engineers. 
Consequently, funding agencies have been reluctant to underwrite the cost of 
practical (large) structural optimization programs, preferring more traditional 
areas of structural mechanics. 

In view of the present situation, the next few years should be dominated 
by increased applications of optimal design, rather than new theoretical develop- 
ments . 


George C. Marshall Space Flight Center 

National Aeronautics and Space Administration - - - - ------ 

Marshall Space Flight Center, Alabama 35812, Aug. 4, 1972 
124-12-11-0000 
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